Strong disorder fixed point in absorbing state phase transitions 
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The effect of quenched disorder on non-equiUbrium phase transitions in the directed percolation 
universality class is studied by a strong disorder renormalization group approach and by density 
matrix renormalization group calculations. We show that for sufficiently strong disorder the critical 
behavior is controlled by a strong disorder fixed point and in one dimension the critical exponents are 
conjectured to be exact: /3 = (3 — \/5)/2 and = 2. For disorder strengths outside the attractive 
region of this fixed point, disorder dependent critical exponents are detected. Existing numerical 
results in two dimensions can be interpreted within a similar scenario. 
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Stochastic many particle systems with a phase tran- 
sition into an absorbing state are of wide interest in 
physics, chemistry and even biology [1]. Much recent 
work has focused on establishing a classification of pos- 
sible universality classes for systems having this type of 
transition [2]. For models with a scalar order parame- 
ter, absence of conservation laws, and short range inter- 
actions, the critical behavior is conjectured to be that of 
directed percolation [3] . Well known models with a phase 
transition in this universality class are the contact pro- 
cess [4] and the Ziff-Gulari-Barshad model of catalytic 
reactions [5]. When there is a conservation law present, 
other universality classes can appear, the best known of 
which is the parity conserving class [6]. 

In this Letter we study the effect of spatially quenched 
disorder on the directed percolation universality class and 
show that for strong enough disorder a new universality 
class of absorbing state phase transitions appears. 

There are two main reasons for performing such a 
study. Firstly, directed percolation has been linked to 
several experimental situations such as catalytic reac- 
tions [5] , depinning transitions [7] , and the flow of granu- 
lar matter [8] (for a review, see [9]). So far however, the 
directed percolation exponents have not shown up in any 
of these experimental realisations. It has been suggested 
[9] that the presence of some form of quenched disorder 
is responsible for this discrepancy. Secondly, while sev- 
eral studies have been performed on the effect of spatially 
quenched disorder [10-15], no coherent picture of these 
effects has been presented. Here, we present such a sce- 
nario by relating the critical behaviour of the strongly 
disordered contact process to that of a class of random 
quantum spin chains, such as the random transverse Ising 
model. Our results mostly concern the one-dimensional 
case, but we will argue that a similar picture also holds 
in d = 2. 

We consider the contact process, in which each site 



of a lattice can either be vacant (0) or be occupied by 
one particle {A). The dynamics of the model is a con- 
tinuous time Markov process in which particles at site k 
can disappear [A — s- 0) with a rate jjLk, while new parti- 
cles can be produced on empty sites (0 — > A), with a rate 
pXk/'i where p is the number of occupied neighbors. Here 
both /ife and \k are independently and identically dis- 
tributed random variables with distributions 7ro(/i) and 
Po(A), respectively. Earlier work on random versions of 
the contact process and of directed percolation focused 
mainly on dynamical aspects, by studying the properties 
of growth starting from a single seed particle. At crit- 
icality, a breakdown of dynamical scaling was observed 
[12], consistent with indications from a field theoretical 
study [13]. The absorbing phase was found to have prop- 
erties similar to that of a Griffiths phase and shows power 
law behaviour with non-universal exponents [12,14]. The 
static critical behaviour of the model has been explored 
less, but exponents in d = 2 have been determined in 
[12]. 

In this Letter we apply for the first time the 'Hamilto- 
nian formalism' [17,18] to the contact process with dis- 
order. This approach has the advantage that it is ideally 
suited to study the asymptotic time regime, which due 
to the slow relaxation properties of the model [15] cannot 
easily be studied with simulations. Since the contact pro- 
cess has the property of duality [19], late and early time 
behaviour can be related, provided suitable initial condi- 
tions are chosen. In this way, the Hamiltonian approach 
can also give insight in the spreading from a seed. 

Consider the contact process on an open chain with L 
sites. The generator of the Markov process is given by 
the 'Hamiltonian' 



Hcp = ^ HkMk + ^ -^{nkRk+i + Rknk+i) (1) 
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in terms of the matrices 
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It is well known [18] that the steady state probabil- 
ity distribution of a stochastic process coincides with 
the ground state of its 'Hamiltonian' and relaxation 
properties can be determined from its low lying spec- 
trum. Although many of the Hamiltonians associated 
with stochastic systems are non-hermitian, techniques 
developed in the study of quantum spin chains can often 
be successfully applied to their stochastic counterparts 
[18,20-22]. 

Here we apply a strong disorder renormalization group 
(RG) approach to the random contact process and the 
results are compared with numerical estimates calcu- 
lated using exact and density matrix renormalization 
(DMRG) based diagonalisations. We first announce our 
basic results which are summarized in the schematic 
phase diagram in Fig. 1 as a function of the control 
parameter A — [lnA]av and the strength of disorder, 
R = [(lnA)^]av — [liiA]av, where, for the sake of simplic- 
ity we put /ifc = 1 (here we use [•]av to denote averaging 
over quenched disorder, while (•) stands for the average in 
the steady state). In the inactive phase, for A < Ac{R) 
the particle density p = limL^oo j Sfc[("'fc)]av is zero, 
whereas in the active phase, for A > Ac{R), p > 0. 
Close to the transition point, for A = — Ac) /Ac i 0, 
the particle density is singular, p ~ A'^ and the cor- 
relation length ^ and relaxation time r diverge as ^ ~ 
[Al^"^^, T ^ |A|~''ii . The anisotropy exponent is defined 
as z = h'\\/i'±. For a semi- infinite system an exponent 
Ps, that describes the singular behavior of the particle 
density at the surface, can be introduced [23] . 

We found that, in one dimension, the universality class 
of the phase transition depends on the strength of disor- 
der, R, and as indicated in Fig. 1 there are two different 
regions of critical behavior. For R > Rc, the critical 
properties of the system are controlled by a strong dis- 
order fixed point, which is situated at i? = oo, and in 
which the critical exponents are conjecturedly exact: 

^(oc)^3_v;5^ /3(-) = l, ,y<f^=2. (2) 

This fixed point is characterised by strongly anisotropic 
scaling 

lnr~^^ (3) 

Thus the anisotropy exponent, is formally infinity. 

For R < Rc, as shown in Fig. 1, there is a line of 
disorder dependent random fixed points (with finite z), 
in which the critical exponents vary continuously with 
R. This behavior is consistent with an extension of the 
Harris criterion for stochastic systems, which predicts 
a cross-over exponent in the non-random system fixed 
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FIG. 1. Schematic phase diagram of the random contact 
process as a function of A and R (see text). The inactive and 
active phases are separated by a phase transition line, along 
which the critical exponents are R dependent for weaker dis- 
order (WD, 7i < Rc), or are determined by a strong disorder 
fixed point at ii = oo, for strong disorder (SD, R > Rc). 

point as ydts{0) = 2/i^_l(0) - 1 > [1,10]. Indeed, ap- 
plying a real space RG method [22] we have calculated 
VdisiO) in good agreement with the Harris criterion. 

In the strong disorder RG, which was introduced [24] 
for random antiferromagnetic spin chains and which was 
later also applied to other one- and higher-dimensional 
systems [16,25,26], one puts the strength of the cou- 
plings and fields in the Hamiltonian (i.e. the rates in the 
stochastic problem) in descending order. The strongest 
one, denoted by fl, sets the energy scale in the system and 
is decimated out and the neighboring rates are replaced 
by (generally) weaker ones, obtained by a perturbation 
calculation. 

In the random contact process, if the fastest rate is 
Xi — fl, then the sites i and i + 1 will either be both 
occupied or both empty almost always, and can thus be 
interpreted as an effective two-site cell state. The ef- 
fective decay rate, p' , to go from AA to 00 can then 
be determined by applying perturbation theory to the 
Hamiltonian (1). At the same time the effective moment 
of the cell state, n', is given by n' = Ui -\- n^+i, where 
in the starting situation Uk = l,Vfc. Conversely, if the 
fastest rate is fXi = fl, site i will almost always be empty, 
and in a perturbative treatment can be eliminated from 
the system. This allows the determination of an effective 
coupling A' between the occupied site i—1 and the empty 
site i + In terms of A = A/VS we obtain the decimation 
equations (k = v^): 

p =K — = , n —Tii + ni+i, A =K . (4) 

Xi Mi 

This decimation procedure is repeated and for strong 
enough disorder, when the probability of generating a 
new rate larger than il is negligible, the energy scale is 
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FIG. 2. Scaling plot of the surface particle density at 
R — 0.5 and at R = 1.5 (inset). 
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FIG. 3. Numerical estimates of the exponents x = /3/i/± 
(circles) and Xa (squares). The broken lines indicate the value 
at the strong disorder fixed point. The errors axe of the same 
order as the size of the symbols. 



continuously lowered and at the same time the proba- 
bihty distributions, 7r(yLt,ri) and P(A, ^2) approach their 
fixed point forms. This RG procedure is easy to im- 
plement numerically. However, in the Id case the fixed 
point distributions can be calculated analytically, using 
the results of a similar analysis for Id random quantum 
spin chains [16,25,26]. Here we mention that according 
to Eq.(4) the renormalised log-/j, rates after large L dec- 
imation steps are expressed as the sum of i + 1 origi- 
nal random log-/i rates minus the sum of L original ran- 
dom log-A rates (plus a negligible constant proportional 
to Iuk). Consequently the fixed point, where the dis- 
tributions 7r(/i, fl) and P{X, ri) are identical, is located 
at [In^Jav = [InAJav Furthermore, at the fixed point 
the central limit theorem implies for the log-energy scale 
— In 17 ~ L^/^, which is equivalent to (3). As a con- 
sequence the probability distributions of the rates are 
broadened without limit during renormalization. Thus 
the ratio of typical neighboring rates goes to either zero 
or infinity. Therefore, in the strong disorder fixed point 
the decimation equations in Eq.(4) become asymptoti- 
cally exact. Thus the critical exponents in this fixed 
point (see (2)), which can be deduced from the analysis 
of the RG equations along similar lines as for the random 
transverse Ising model [16], are presumably exact, too. 

In the following we check these predictions by numer- 
ical calculations on lattices with L up to 32 with the 
DMRG method. For the disorder we used a bimodal dis- 
tribution with — I and P{X) — [S{X — A+) -f- ^(A — 
A_)]/2, where A± — exp{A ± VR), so that A and R are 
the parameters used in Fig. 1. For smaller sizes, L < 14, 
we performed exact disorder averages, whereas for larger 
L-s we considered at least 10000 samples for each size. 

In the actual calculation we set fJ.L = and deter- 
mine the average particle density profile, pk = [{nk)]a.v, 
so that the order-parameter in the bulk (at the sur- 
face) is defined as p = pL/2 {ps = Pi)- Close to the 



critical point, the densities obey scaling laws, such as 
p(L,A) = L-^p(ALi/''^), where x = A sim- 

ilar relation holds for pa, with the surface exponent 
Xg — Ps/v^. From an optimal data collapse, as illus- 
trated in Fig. 2, the position of the transition point, Ac^ 
and the critical exponents are calculated. 

For varying disorder strength, i?, the estimated crit- 
ical exponents, x and plotted in Fig. 3, where 
one can observe the two regions announced in Fig. 1. 
The crossover occurs near Rc ~ 1.5. For higher disorder 
strenghts the critical exponents seem to saturate at the 
values given in (2). The estimate for v±^ = 1.7 ± .3 at 
R = 1.5 is also consistent with v^^\ The error on this 
exponent is larger because of stronger finite size correc- 
tions. 

We note that the same type of weak-to-strong disor- 
der cross-over scenario has been observed in a class of 
random quantum spin chains [27]. For these models, the 
strong disorder RG analysis leads to decimation equa- 
tions of the same form as Eq.(4), with also k > 1. For 
these K-values, some of the generated couplings (rates) 
can be larger than the decimated one. For weak enough 
disorder, this happens so often that the RG approach 
does not work. For stronger, but still finite disorder, 
however, such type of decimation steps are rare so that 
during renormalization the system flows into the strong 
disorder fixed point [27]. 

We have also studied the dynamical scaling of the 
model by calculating the probability distribution of the 
gap F between the groundstate and the first excited state 
of Hcp- In the small F limit this distribution behaves 
as P(lnr) ^ F^/^ . Here z' is the disorder induced 
dynamical exponent, and the true dynamical exponent 
of the system is given by z = max(z', z(0)). In the 
strong disorder fixed point z' is formally infinity and 
the scaling behavior of the gap distribution is given by 
P(lnr) = p(\nT/VL) /Vl. As shown in Fig. 4 the 
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FIG. 4. Scaling of the gap to the first excited state at 
R = 1.5 and A = Ac = 1.18 {L = 16 (-I-), 20 (□), 24 Q)- 



above scaling form is indeed satisfied in the strong disor- 
der region at R= 1.5. 

Based on the decimation rules in Eq.(4), which are very 
similar to that of the random transverse Ising model, 
we expect that even in higher dimensions these ran- 
dom stochastic models will renormalize into the strong 
disorder fixed point for the higher dimensional random 
transverse Ising model [28]. Indeed from [28] we obtain 
uj_ = 1.07±.15 and /3 = 1.07±.12, consistent with the val- 
ues obtained for a contact process with sufficiently strong 
dilution [12]. For smaller values of the dilution, /3 seems 
to vary continuously [12], indicating that a scenario sim- 
ilar to that in (i = 1 is valid. Moreover, the breakdown 
of dynamical scaling observed in [12] can be related to 
the strongly anisotropic scaling (3). Using tp = A2± .06 
[28] , it is even possible to predict [29] the values of some 
of the dynamic exponents introduced in [12]. 

In conclusion, we have found that the critical proper- 
ties of absorbing states models in the directed percolation 
universality class in the presence of random transition 
rates are controlled by a strong disorder fixed point, if 
the disorder is sufficiently strong. The exponents associ- 
ated with this transition, are conjectured exact. This is 
the first example of an absorbing state phase transition 
with quenched disorder for which exact information be- 
comes available. A more extensive paper on our results 
will be published elsewhere [29]. 
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